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The statistically steady humidity distribution resulting from an 
interaction of advection, modeled as an uncorrelated random 
walk of moist parcels on an isentropic surface, and a vapour 
sink, modeled as immediate condensation whenever the specific 
humidity exceeds a specified saturation humidity, is explored with 
theory and simulation. A source supplies moisture at the deep- 
tropical southern boundary of the domain, and the saturation 
humidity is specified as a monotonically decreasing function 
of distance from the boundary. The boundary source balances 
the interior condensation sink, so that a stationary spatially 
inhomogeneous humidity distribution emerges. An exact solution 
of the Fokker-Planck equation delivers a simple expression for 
the resulting probability density function (PDF) of the water 
vapour field and also of the relative humidity. This solution 
agrees completely with a numerical simulation of the process, 
and the humidity PDF exhibits several features of interest, such 
as bimodality close to the source and unimodality further from 
the source. The PDFs of specific and relative humidity are broad 
and non-Gaussian. The domain averaged relative humidity PDF 
is bimodal with distinct moist and dry peaks, a feature which we 
show agrees with middleworld isentropic PDFs derived from the 
ERA interim dataset. 

Copyright (c) 2010 Royal Meteorological Society 

Key Words: advection, diffusion, condensation, probability density function, relative humidity, 
specific humidity, water vapor, bimodal 



Received . . . 
Citation: . . . 



1. Introduction 

Water vapour plays a significant role in diverse 
problems pertaining to the dynamics of the Earth's 
climate (Pierrehumbert 2002; Schneider, O'Gorman 
and Levine 2010). In particular, with regards to 
its radiative effects, estimating the amount of water 
vapour in the upper troposphere is crucial (Held and 
Soden 2000; Spencer and Braswell 1997). Significant 
efforts — especially directed towards the subtropical 
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troposphere — have yielded a key insight into 
the nonlocal nature of processes that control the 
distribution of water vapour (Yang and Pierrehumbert 
1994; Sherwood 1996; Salathe and Hartman 1997; 
Pierrehumbert 1998; Pierrehumbert and Roca 1998; 
Galewsky, Sobel and Held 2005; Brogniez, Roca and 
Picon 2009). 

The resultant framework for understanding the 
evolution of water vapour, or generally any condensable 



2 



J. Sukhatme and W. R. Young 



substance in a fluid dynamical setting, is known as the 
advection-condensation (AC) model. The recent review 
by Sherwood et. al (2009) provides an introduction, 
while the book chapter by Pierrehumbert, Brogniez 
and Roca (2005) (hereafter, PBR) provides a thorough 
overview of the AC model. 

Further, as elaborated on by Spencer and Braswell 
(1997) and PBR, in addition to the amount, the 
probability density function (PDF) of the water vapour 
field has a strong influence on the longwave cooling 
of the atmosphere. Explicit calculations illustrating 
the effect of the PDF on the outgoing longwave 
radiation (OLR) are given by Zhang, Soden and Mapes 
(2003). The importance of the PDF stems from the 
approximately logarithmic dependence of the change 
in OLR to the specific humidity of the water vapour 
in the domain. Further discussion of this logarithmic 
dependence which results from a spectral broadening 
of absorption peaks can be found in the text by 
Pierrehumbert (2010). 

Following Soden and Bretherton (1993), there 
have been numerous efforts that examine water vapour 
PDFs in the tropical and subtropical troposphere. A 
principal feature noted by these studies is the non- 
normality of the PDFs. A range of distributions, 
from bimodality in the deep tropics (Brown and 
Zhang 1997; Zhang, Soden and Mapes 2003; Mote 
and Frey 2006; Luo, Kley and Johnson 2007), to 
unimodal PDFs with a roughly lognormal or power-law 
form in the subtropics have been documented (Soden 
and Bretherton 1993; Sherwood, Kursinski and Read 
2006). Ryoo, Igusa and Waugh (2009) have emphasized 
that a salient and baffling feature of water vapour 
PDFs is their spatial inhomogeneity. Models based 
on the statistics of subsidence drying and random re- 
moistening of parcels reproduce many of the observed 
features of the humidity PDFs, provided that one 
considers the model parameters to be functions of 
position determined to match observations (Sherwood, 
Kursinski and Read 2006; Ryoo, Igusa and Waugh 
2009). 

One of the main achievements of this paper is 
an exact solution for the water vapour PDF that 
characterizes the statistically steady state emerging 
from an interaction of advection, condensation and 
a spatially localized moisture source. In particular, 
we consider moist parcels advected by a white-noise 
stochastic velocity field. Each parcel carries a specific 
humidity and experiences condensation whenever and 
wherever the specific humidity exceeds a specified 
saturation profile. The parcel humidity is reset to 
saturation by a moist source at the southern boundary 
of the domain. In the ultimate statistically steady state, 
the southern vapour source is balanced by condensation 
throughout the interior of the domain. In contrast to 
earlier considerations of the AC model with white-noise 
advection (O'Gorman and Schneider 2006), we consider 
statistically steady states, rather than the initial value 
problem, and we obtain the full water- vapour PDF 
rather than the mean humidity. 

The AC model is described in more detail in 
sections 2 and 3, and in section 4 we solve the Fokker- 
Planck equation to obtain the PDF of the specific 
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humidity. The utility of the PDF is demonstrated by 
calculating the average moisture flux and condensation. 
In section 5 we focus on a particular saturation 
profile — the standard exponential model — and 
exhibit analytical solutions for the specific and relative 
humidity PDFs; these are shown to be in very good 
agreement with PDFs estimated from a Monte Carlo 
simulation. In section 6 we analyse daily isentropic 
relative humidity data from the ERA interim product, 
and show that the PDFs derived from this dataset 
agree with principal features predicted by our idealized 
model. A discussion of the PDFs, with implications for 
the atmospheric radiation budget and the present-day 
distribution of water vapour, followed by a conclusion 
section 7 ends the paper. 

2. The AC model 

The central quantity in the AC model is the specific 
humidity, defined by 

def mass of water vapour in a parcel 
' total mass of the parcel ' 

where x = {x^y) is position on an isentropic surface. 
Our main focus here is the extratropical vapour 
distribution, which is determined in large part by the 
wandering of air parcels on isentropic surfaces and 
the temperature changes associated with latitudinal 
excursions. 

Following PBR, the AC model is: 

dtq^U'Vq = S -C. (2) 

Above, S is the vapour source, C is the condensation 
sink and u{x^t) is velocity on an isentropic surface. In 
the absence of sources and sinks, Dq/Dt = 0, i.e. the 
specific humidity is a materially conserved quantity. 
The AC formalism does not consider the diffusive 
homogenization of water vapour e.g., via the addition 
of x^^Q oil the right of (2). This limitation of the model 
is discussed further in section 6. 

Latent heat release and radiative cooling result in 
significant cross-isentropic motion in the troposphere, 
and we have ignored this process in formulating a two- 
dimensional isentropic model in (2). We emphasize that 
the AC formulation is not inherently limited to two- 
dimensions or to motion on isentropic surfaces. But for 
ease of exposition, and simplicity, isentropic motion is 
a useful starting point. 

The condensation sink on the right of (2) might 
be modeled with 

C = [q{x, t) - q,{x)] H [q{x, t) - q,{x)] , (3) 

where Tc is the timescale associated with condensation, 
qs{x) is the prescribed saturation specific humidity; 
H{x) denotes the Heaviside step function, which is zero 
if the X < and one if x > 0. Observed atmospheric 
supersaturations rarely exceed a few percent and 
therefore the sink C is extremely effective (Wallace 
and Hobbs 1977). In other words, we are concerned 
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with the rapid-condensation hmit Tc ^ 0. In the rapid- 
condensation hmit, the details in C are unimportant 
and condensation can instead be represented as a rule 
enforcing subsaturation: 



q{x, t) min [q{x, t), q^{x) 



(4) 



Therefore, qix^t) < q^ix) and there is an upper bound 
on the specific humidity at every x. Globally, the 
maximum tolerable specific humidity is 



def 



max(7s(ic) • 

Va3 



(5) 



If the fiow is ergodic then a wandering parcel 
explores the entire domain and eventually encounters 
the dryest regions, with the smallest value of qs{x). 
Thus in the absence of sources (Sukhatme and 
Pierrehumbert 2006) 



lim q{x^ t) 

t^oo 



where 



def . / X 

= mmqs[x) . 



(6) 



(7) 



On the other hand, without advection, the source will 
eventually saturate qix^t) at every point x where S ^ 
0. Therefore, interesting statistically steady states are 
expected only if there is both a source S and advection 
n(ic, t). 

The problem is to determine the evolution and 
the statistics of the specific humidity q{x,t)^ given the 
velocity field u{x^t) and the saturation profile qs{x). 
Due to its physical significance, the principal statistical 
feature we focus here on is the probability density 
function (PDF) of q. Quantities such as the mean 
moisture fiux and the mean condensation rate are also 
crucial, and have been examined in the AC framework 
(O'Gorman and Schneider 2006; O'Gorman, Lamquin 
and Schneider 2011). These mean quantities can be 
deduced from PDF of q via integration over q. 



3. Moist brownian parcels 

The simplest representation of random transport is 
obtained by taking the velocity u{x^ t) in (2) as having 
rapid temporal decorrelation so that the motion of an 
air parcel is brownian. Then the spatial evolution of an 
ensemble of parcels, each carrying its own individual 
value of g, is described by the diffusion equation. We 
also assume that depends only on the meridional 
coordinate y in the domain < y < L, and that q^ 
decreases monotonically from ^max at y = to ^min 
at y = L. Each parcel is located by a point in the 
configuration space {q^y). Because of the condensation 
rule (4), the ensemble of parcels is confined to the 
shaded region of configuration space illustrated in 
Figure 1. 

The AC model is then expressed as a set of coupled 
stochastic differential equations governing the position 
Y{t) and specific humidity Q{t) of a moist brownian 
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Figure 1. In the rapid-condensation limit the accessible domain 
in the (q, y) plane is the shaded region shown above. The model 
assumes that Qsiy^ decreases monotonically from Q'max — Qs 

(0) to 

Qmin = Qs(L). The functions qs(y) and ys(q) are inverses of each 
other: qs{ys{q)) = q. 



parcel 

dY{t) = V^dW(t), 

dQ{t) = {S[Y{t)]-C[Y{t),Qm dt (8) 

where W{t) is a Weiner process, is the brownian 
diffusivity associated with the random walk, and S 
and C are the source and condensation functions in (2). 
The numerical simulations described below amount to 
Euler-Maruyama discretization of (8). The midlatitude 
baroclinic eddies responsible for transport have non- 
trivial spatio-temporal correlations (Sukhatme 2005; 
O'Gorman and Schneider 2006). Thus the white noise 
advection used in (8) and in earlier AC studies is an 
idealization. 

The statistical properties of the ensemble of moist 
brownian parcels are characterized by a PDF P{q^ y^ t), 
such that the expected number of parcels in dqdy is 
equal to 

NP{q,y,t)dqdy, (9) 

where N is total number of parcels. 

The evolution of P{q,y^t) is governed by the 
Fokker-Planck equation 

dtP + dg[{s - C)P] = KbdlP , (10) 

with the normalization. 



P{q,y,t) dqdy = l. 



(11) 



In the limit of rapid condensation, as gmin ^ 
Q{y^t) < qs{y), the normalization takes the form 

/ / P{q,y,t)dqdy = l. (12) 



Integrating P{q^y,t) over q produces the marginal 
density 

(13) 

Qmin 



p{y,t) = / P{q,y,t)dq. 
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Because the brownian motion of a parcel is independent 
of the specific humidity it carries, we expect that 
integration of (10) over q wih show that p{y^t) 
satisfies the diffusion equation, and indeed this result 
is immediate 



dtP = i^bdyP . 



(14) 



If parcels are refiected back into the domain at y = 
and L then the solution of (14) is limt^oo Piv^t) = L~^. 
Thus if we seek long-time equilibrium solutions of the 
Fokker-Planck equation (10), then we require 



L- 



(15) 



Note that the left of (15) is consistent with the global 
normalization in (11). In the limit Tc ^ 0, the steady, 
marginal normalization condition (15) applies at every 
7/, and for all models of S and C. In other words, in a 
statistical steady state the normalization condition (15) 
is equivalent to enforcing the rapid-condensation rule 
in (5). The normalization (15) figures prominently in 
the following solution of the Fokker-Planck equation. 
Similar "variable limit" normalization conditions are 
encountered in the solution of integrate-and-fire neuron 
models (Fusi and Mattia 1998). 

4. Steady solution with "resetting" at y = 

Now we adopt a "resetting" model for the vapour 
source S: on encountering the southern boundary at 
7/ = 0, parcels are refiected back into the domain with 
a new value of q chosen from a specified probability 
density function ^(g), with the normalization 



Q'max 

^q)dq = l. 



(16) 



For example, completely saturating the parcels at 
y = corresponds to the choice ^{q) = S{q — ^max)- 
This is a simple representation of tropical moistening. 
Resetting boundary conditions have been employed in 
a previous study of AC on isentropic surfaces (Yang 
and Pierrehumbert 1994), as well as in advection- 
diffusion of passive scalars (Pierrehumbert 1994; 
Neufeld, Haynes and Picard 2000). 

We emphasize that resetting the humidity only at 
7/ = is a strong model assumption that distinguishes 
this work from the AC models developed by Sherwood, 
Kursinski and Read (2006) and Ryoo, Igusa and Waugh 
(2009). The statistical model developed in those works 
assumes that random resetting occurs throughout the 
domain as an exponential or gamma stochastic process 
with a characteristic time scale rMoist(^c)- In effect, 
by resetting only at ^ = 0, we are considering the 
extreme case in which rMoist(^^) = co, except at = 
0. Physically, the picture we have in mind is that 
parcels experience moistening when they encounter 
convective regions which are restricted to the tropics, 
and therefore y is the distance from this saturated zone. 
In fact, this serves as another limiting case for the AC 
model along with the re-setting throughout the domain 
scenario considered by Sherwood, Kursinski and Read 
(2006) and Ryoo, Igusa and Waugh (2009). 



An important consequence of resetting only at 
y = is that the dryest parcels, which are created at 
y = L with q = g'min, maintain their extreme dryness 
till they eventually strike y = 0. The solution below 
requires explicit consideration of these dryest parcels 
by inclusion of a component S{q — ^min), referred to as 
the "dry-spike", in the PDF. 

4.1. Solution of the Fokker-Planck equation 

In the rapid-condensation limit, the accessible part of 
the (g, ?/) -space is the shaded area in Figure 1, where 
C = 0. But with the resetting boundary condition, the 
source S is also zero within the shaded region. Thus 
the steady Fokker-Planck equation in the shaded region 
of Figure 1 collapses to = 0, with the immediate 



solution 



P{q.y) = ^^yB{q). 



(17) 



The first term on the right of (17) is determined 
by the resetting condition and by application of the 
normalization in (15) at = 0. 

Corresponding to the dryest parcels in the domain, 
B{q) in (17) must contain a component proportional to 
^^(^ — ^min)- This "dry spike" contains parcels created 
by encounters with the northern boundary at y = L. 
These considerations refine (17) to 



P{q,y) 



L 



Sjq- gmin) 
L2 



F{q) 



y- (18) 



The term L '^S{q — qmin) in (18) is obtained by 
applying the normalization condition (15) ai y = L. 

In (18), F{q) is the smooth part of the PDF, which 
is determined by substituting into the normalization 
(15). This requirement leads to 

- ^q')dq' = l + L^ F{q')dq'. (19) 

Inspecting (19), we see that it is possible to solve for 
F{q) by using q^ rather than y^ as the independent 
variable. That is, following the notation in Figure 1, 
we use y = ys{q) to re- write (19) as 



L 

ys{q) Jq 



/ <^>{q')dq' = l + L^ / F{q')dq'. 



(20) 



The derivative with respect to q then determines F{q) 
i.e. 

d A{q) 



LF{q) , , , , 
dqvM) 

where A(g) is the cumulative distribution 



def 



$(q')d</' 



(21) 



(22) 



Assembling the pieces, we obtain the PDF 



1 dA 

'L~dq 



y_ 

L 



S{q 



d A{q) 



dq ysiq) 



(23) 
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The expression above apphes only within the 
shaded region in Figure 1 where q < ^s(^); in the 
supersaturated region q > qs{y)^ P{q, y) = 0. The PDF 
is discontinuous at the saturation boundary q = qs{y)- 
This completes the solution for the steady-state Fokker- 
Planck equation with resetting forcing at ^ = 0. 

In summary, in the rapid-condensation limit, 
and with brownian motion, the steady Fokker-Planck 
equation is solved exactly by (23). The normalization 
constraint in (15), and the use q as the independent 
variable, produces this relatively simple solution. 
In fact, the normalization constraint (15) can also 
be regarded as a differential boundary condition. 
Manipulating this boundary condition provides an 
alternate route to the solution outlined in Appendix 
A. 

4.2. The global PDF 



4.4' Average humidity and condensation rate 

To understand the transport of moisture in this model, 
and the associated mean condensation rate, we return 
to (10) and momentarily retreat from the rapid- 
condensation limit*. Then, multiplying (10) by q and 
integrating from q^nin to ^max, we have 



qdq(CP) dq = 



d?/2 



(30) 



where the condensation sink C is given by the model 
in (3). Using this expression for C and integrating the 
left of (30) by parts, we obtain the mean condensation 
rate, defined by 



{q - ^s) P{q^ y) M , 



(31) 



The global PDF of specific humidity is the marginal 
density 

(24) 



9{q) = / P{Q^y)(^y- 

Jo 



The integral on the right can be evaluated using the 
solution in (23), and one finds 

9{q) = ls{q - g„i„) -~ {ys{q)Hq)) ■ (25) 

A comforting check on (25) is that g{q) satisfies the 
normalization 



g{q) dq = l, 



(26) 



which is consistent with (12). 

Independent of model details, such as the 
specification of qs{y) and A(g), the dry spike, 6{q — 
Qmin) in (25), contains exactly half of the parcels 
in the ensemble. Prof. P. O'Gorman noted that this 
feature can be understood on the basis of symmetry: 
exactly half the parcels have more recently encountered 
the moist southern boundary than the dry northern 
boundary. The former half have q > g'min, while the 
other half constitute the dry spike with q = qmin- 

4.3. Averages 

Knowledge of the stationary PDF allows us to evaluate 
the average of any function of say /(g), as 



def 



{f){y) = L f{q)P{q,y)dq. 



(27) 



dy 



r 



d/ A(g) 
dq ys{q) 



dq. 



(29) 



(C) = Kb 



dHq) 
dy2 



(32) 



The right hand side of (32) is the convergence of the 
moisture flux, which balances the mean condensation 
(C) on the left. 

This argument shows that the flux of moisture has 
the same diffusivity Kb as that of a brownian parcel. 
As explained by O'Gorman and Schneider (2006), this 
is a special property of the brownian limit: if the 
velocity has a nonzero integral time scale (e.g., as in 
the Ornstein-Ulenbeck process) then the diffusivity of 
moisture will be systematically less than the random- 
walk diffusivity. We emphasize that while the moisture 
flux is diffusive, even with a Brownian advecting 
flow, the condensation sink cannot be adequately 
represented by a bulk diffusivity. This inadequacy of 
a mean-field representation of condensation is clearly 
demonstrated in PBR. 

In the rapid condensation limit, the right of (32) 
can be computed by putting f{q) = q in (29). Thus if 
Tc ^ 0, the mean condensation rate is 



Kb dQs 

y dy 



A(9s(y)) . 



(33) 



Using the solution for P{q^y) in (23), one can show 
that the average defined above is 

dA df A(q) 

(/) = /(^s)A(^s) - / m-^dq-y / /-p(d^. 

(28) 

Taking the derivative of (28) with respect to y we 
obtain, after some remarkable cancellations. 



Although C (X ^ 00, the mean condensation (C) in 
(33) is independent of Tc. This singular limit is achieved 
because the nonzero probability of supersaturation is 
confined to a boundary layer, lying just above the 
saturation boundary q = qs{y) in Figure 1. 

If the resetting PDF ^{q) is non-singular at ^max, 
then (C) in (33) is finite as ^ 0, despite the factor 
y~^. Specifically, one can show that liiRy^o A{qs) / y = 
-^(gmax)(d(7s/d^). 

5. An example 

The solution from the previous section is most simply 
illustrated by supposing that the resetting produces 



*We continue using resetting forcing, so that *S = when y > 0. 
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complete saturation at 7/ = 0, 



(34) 



so that A{q) = 1. For the saturation humidity we use 
the exponential model 

=9maxe-"^. (35) 

This form of has been employed previously in 
AC studies (Pierrehumbert, Brogniez and Roca 2005; 
O'Gorman and Schneider 2006). 

With ^{q) and q^iy) in (34) and (35), the steady 
solution of the Fokker-Planck equation is 

^{q- ^max) 



and away from the southern boundary 



(36) 



Piq,0<y) = 



Siq 







a 



L 



(7 log (g/gmax) 



(37) 



The solution (37) can be non-dimensionalized with 

def 



def r 

= aL . 



def / 
g* — Q/ ^max 7 



y, = y/L, (38) 

and P*(<7*,^*) = qniaxLP{q^y). This scaling identifies 
the fundamental non-dimensional control parameter 
a*, and the scaled PDF is 



P*(g*,0 < y^) = y^ 



S{q^ - e) + 



g^log (g*) 



with 



def <7n 



(40) 



At the southern boundary P^(g^,0) = S{q^ — 1). The 
main point is that the shape of P^{q^,y^) is controlled 
by a single parameter a^. 

There is an interesting change in the structure of 
at a* =2: is bimodal for all y^ifa^ < 2; if > 2 
then P* is bimodal only in the region < y^ < 
This is illustrated in Figure 2 which shows the smooth 
portion of P* with a* = 4 and = 1. The PDF with 
a* = 4 changes from a bimodal to a unimodal form as 
one crosses y^ = = 0.5. Indeed, for y^ > 0.5 the 
a* = 4 PDF falls monotonically from q^ = e. In the 
case with a* = 1, the PDF increases with increasing 
q^ at all values of Therefore, in combination with 
S{q^ — e), the PDF is always bimodal if a* < 2. 

On a middleworld isentrope, the saturation mixing 
ratio usually changes by more than an order of 
magnitude as one spans the entire midlatitudes, hence 
from (35), we expect that a* > 2 is the relevant 
atmospheric case. Therefore, close to the vapour source 
the PDF is bimodal. One mode, the dry spike 6{q — 
gmin), corresponds to a singularity at gmin, after which 
the PDF falls off, only to rise again to a secondary 
moist peak as g ^ Qsiy)- Further from the vapour 
source (i.e. for y^ > the PDF is unimodal with 

a dry spike at ^min followed by a gradual monotonic 
decay to a minimum at g = qsiy)- 

We proceed to discuss the implications of the 
solution in (37) using dimensional variables. The non- 
dimensional versions of all formulas are obtained with 
Q, — y Q^>K and gmin — ^ 




Figure 2. The behaviour of the smooth portion of P*(g*,y*) 
for two choices of q;>k 5 we do not show the dry spike at Q'min- 

The 

major axes correspond to o;* = 4; there is a transition from a 
bimodal to unimodal PDF at = = 0.5. Larger ?/*'s have 

smaller ranges of permissible , and hence the curves are ordered 
such that increases as we move to the left. The transition 
curve at = 2/q;* = 0.5 is shown by the dark line in the plot. 
The inset corresponds to o;* = 1. The curves are ordered as in 
the main plot. In this case, counting the dry spike at gmin, the 
PDF is bimodal at all values of y*. 



(39) 5.1. Statistics in a strip yi < y < y2 



The expression in (36) and (37) provides the PDF at a 
specified y. From an observational and numerical view, 
recording the PDF at a particular y can be difficult. It 
is more likely that one obtains an estimate of the PDF 
over a strip, i.e. for y G [^1,^2]- Integrating (37) over 
such a strip, and denoting this by Pi2(^), we have 



ry2 

Pi2{q) = / P{q,y)dy, 
' S{q-qmin) 

where the strip function is 



2L^log ((^/^max) 



(41) 



[yiiQ) 



if Qmin <q< qs{y2); 
if qs{y2) <q< qs{yi)- 



(42) 



The slightly complicated form of Si2{q) is because when 
qs{y2) < q < qs{yi)^ the upper limit of the integral is 
Vsiq) < y2' We compare P12 above with the results of 
a Monte Carlo simulation in Figure 3. 

The effect of the strip function 512 is to clip 
the secondary peak at q = qs^ and this clipping 
becomes more extreme if 7/2 — yi is made larger. Thus, 
aggregating measurements of q from different spatial 
locations will obscure the bimodal structure of P{q^ y). 

5.2. The global PDF of q 

In the extreme case ^1 -> and y2 ^ we obtain the 
global PDF of the specific humidity, denoted by g{q)^ 
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Figure 3. Comparison between Pi2(q,y) from Monte Carlo 
with the smooth portion of the analytic result in (41). In this 
illustration a* = 5. Crosses, open circles and stars are the Monte- 
Carlo PDFs while the smooth solid curves are the analytical 
expression. Top panel: The bimodal PDFs for < 2/q;*. Middle 
panel: The PDFs for larger which show the transition from a 
bimodal to a unimodal distribution. Bottom panel: The unimodal 
PDFs far from the source, i.e. for > 2/q;*. 



in (24). Using (25), the global PDF in this example is 
1 .1 



9{q) = 2^{Q 



+ 



2aLq ' 
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(43) 



where q G [qmin 
PDF is g{q) - q-^ . 



The smooth part of the global 



5.3. Monte Carlo Simulations 

To test the analytic solution above, we proceed to a 
Monte Carlo simulation of our system. Typically our 
simulations use N = 10^ parcels which are released in 
the domain y G [0, L]. We have used a = 1, L = 5 which 
gives a* = 5. This implies a change of approximately 
two orders of magnitude in the saturation specific 
humidity across the domain. The ^/-coordinate of the 
parcels are initially uniformly distributed and are 
tagged with the minimum saturated specific humidity, 
i.e. on parcel z, qi{0) = <7min- The parcels perform an 
uncorrelated random walk in and parcel i refiected 
back into the domain if yi > L or if yi < 0. The 
forcing S is implemented by setting qi = g'max when 
a parcel is refiected at y = 0. Rapid condensation is 
implemented by enforcing q^ = mm[qi, QsiVi)] after each 
random displacement. 

To approximate brownian motion, the root mean 
square step length of the random walk, should be as 
small as computationally feasible — we typically used 
£ = 2 X 10~^L. We specify the brownian diffusivity 
as hib = 0.5, and then the time between steps, r, is 
computed via 

£^ 



2r 



(44) 



The simulation runs till P{q, y) attains a stationary 
form and then data is collected in several strips for 
comparison with (41). 

The results of the simulations are shown in 
Figure 3. The top panel of this figure shows a 
comparison between the smooth portion of the 
theoretical expression in (41) and the numerically 
estimated PDF for < 2/a*. There is good agreement 
between the Monte- Carlo simulation and the Fokker- 
Planck solution. Some minor discrepancies are evident 
as g ^ ^min; this is anticipated as we have suppressed 
the numerical spike q = q^^^ (indeed, the theoretical 
PDF is given by a ^-function at this location). The 
middle panel of Figure 3 shows the PDF at larger 
7/, specifically around the region where the PDF 
transitions from a bimodal to a unimodal form. The 
bottom panel of Figure 3 shows the PDF's far from the 
source, where the PDFs have a single peak q = q^im 
followed by decay and termination dXq = qg (y) . Finally, 
the first panel of Figure 4 shows the global PDF g{q) 
in (43). As per the log-log inset, the global PDF closely 
follows the result in (43), with a spike dit q = gmin 
followed by an algebraic tail g{q) ^ q~^ • 

5.4- Transport and condensation 

To understand the water vapour distribution, it is 
important to quantify the transport of moisture. In the 
context of the AC model, it was shown by O' Gorman 
and Schneider (2006) that brownian motion results in 
the mean moisture fiux 
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dy 
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Figure 4. Top panel: The global PDF of the specific humidity 
g(q) from Monte Carlo simulations. The inset shows a log-log 
plot which reveals a spike at gmin followed by an algebraic 
tail as predicted by (43). Bottom panel: The analytical curve h{r) 
from (51) (solid line) and the global RH PDF from the numerical 
simulation (dotted circles). 



where {q) is the mean specific humidity and /^^ is the 
brownian diffusivity. Then, as in (32), in a statistical 
steady state the convergence of (J^) balances the mean 
condensation (C) at every y. 

Using f{q) = q in (29), and qs{y) in (35), the 
average gradient of specific humidity is obtained as 



dy 



where E{x) =^ J^{e~'^/t) dt is the exponential integral. 
The moisture flux {J^) can also be diagnosed from 
the simulation and thus (45) is tested with agreement 
shown in Figure 5. 

5.5. The Relative Humidity (RH) 

It is the immense subsaturation of the troposphere 
that makes the study of water vapour an interesting 
and challenging problem (Spencer and Braswell 



[E(ai) - E(ay)] , 



(46) 



Figure 5. A comparison of the flux {J^) measured in the Monte 
Carlo simulation with the theoretical result using (45) and (46). 



1997; Pierrehumbert, Brogniez and Roca 2005). This 
subsaturation is most starkly revealed by examining 
the RH, which is 

def q ..^x 

Qs{y) 

Converting the PDF in (37) to a PDF for r, denoted 
by Q(r, y), gives (for ^ > 0) 



Q{r^y) = 



S{r 







a 



, (48) 



L rlog^(er/rmin). 
where r G [rmin(y), 1], e is defined (40), and 

rminiy) = ee"« . (49) 

A comparison between (48) and the Monte Carlo 
simulation is shown in Figure 6. The main mismatch we 
see is for r ^ 1, where the numerical PDF's taper off 
while expression (48) continues to rise till r = 1. This 
discrepancy is due to the collection of numerical data 
in a strip (as before), while (48) provides an expression 
for the PDF at a particular location. Note that the 
numerical PDF's show a spike at Vminiy) — shown 
only in the first panel — which agrees with the shifting 
^-function in (48). Bimodality of the PDF at small 
y is evident, while farther from the source there is a 
peak at rmin(^) (suppressed in the plot), after which 
the PDF gradually levels off as r ^ 1. The transition 
from a bimodal to a unimodal PDF, as well as the 
PDFs far from the source, are shown in the second 
and third panels of Figure 6 respectively. 

5.6. The global RH PDF 



The global PDF, h{r), of RH is 



h{r) / 
Jo 



log(r/e)/a 

Q{r, y)dy. 



Evaluating (50) using Q{r^y) in (48) gives 



h{r) = 



1 



aLr 



log 



aL 



log(l/r) 



(50) 



(51) 
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As a check one can verify that on integration over the 
range of RH, that is e < r < 1, the global PDF in (51) 
normalizes to unity. The lower panel of Figure 4 shows 
a comparison between the Monte Carlo global RH PDF 
and (51). There is a good match between the theoretical 
prediction (50) and the simulation. The global PDF 
exhibits a mix of the distinct features of the PDFs at 
different y e.g., h{r) is bimodal with a primary peak 
at small RH, followed by a decay at intermediate RH, 
and then a secondary "moist peak" as r — > 1. 

6. Discussion 

6.1. Radiative effects of water vapour 

In the introduction, we referred to the importance 
of the vapour PDF in understanding the effect of 
humidity on the OLR. In fact, the effect of humidity 
on OLR was quantified by Zhang, Soden and Mapes 
(2003) in terms of the PDF of the RH field. 
Their numerical experiments, employing a full-fledged 
radiation code, showed that a domain with a bimodal 
RH PDF allowed for significantly larger cooling than 
one with a unimodal distribution (holding the total 
and mean the same in the two cases). In the present 
context, we see that the PDF of the RH over the 
domain, i.e. the lower panel of Figure 4, is bimodal. 
Therefore, along with regions of moderate RH, the 
domain consists of subsets with very small RH (the first 
peak in the aforementioned plot) and those that are 
nearly saturated (the second peak). A clearer picture 
is seen in the upper panel of Figure 7 which shows 
the average specific humidity and RH over the domain, 
i.e. {q){y) and {r){y) respectively. While the specific 
humidity decays monotonically from the source, the RH 
has a pronounced mid-domain minimum representing 
a pool of "relatively" dry air which is flanked on either 
side by nearly saturated regions. 

As described in PBR, the presence of water vapour 
affects the OLR in an approximately logarithmic 
fashion, i.e. OLR ex — log(g'). Therefore, the importance 
of knowing the PDF of q at every y is illustrated 
by comparing log((g)) with (log(g)) — it is the latter 
quantity that accounts for fluctuations in the specific 
humidity at every location and requires a knowledge of 
P{q^ y). As expected (see for example the discussion in 
PBR), and as is shown in the lower panel of Figure 
7^ (especially the inset, which shows the difference 
between the two estimates) the fluctuations in the 
specific humidity field increase the OLR as compared 
to the estimate based only on the mean profile {q){y). 



Figure 6. Comparison between Qi2(T,y) from Monte Carlo 
with the smooth portion of the analytic result in (48). In this 
illustration q;^ = 5 and 1/2 — yi = 0.02L. Crosses, open circles 
and stars are the Monte-Carlo PDFs while the smooth solid 
curves are the analytical expression. Top panel: The bimodal 
PDFs for < 2/q;* (the curve in the middle corresponds to 
= 0.08). Middle panel: The PDFs for larger which show the 
transition from a bimodal to a unimodal distribution. Bottom 
panel: The unimodal PDFs far from the source, i.e. for > 



6.2. Bimodality of the RH PDF 

It is interesting that the global PDF of q in the upper 
panel of Figure 4 is unimodal, while the global PDF 
of RH in the lower panel of Figure 4 is bimodal. Thus, 
according to our model, RH bimodality is so robust 



"'"Both curves start at zero as we have chosen qn 
the value at y = is log(gmax)- 



1, in general 
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Figure 7. Upper panel: The average profiles of the specific 
humidity and RH computed from the PDFs in (37) and (48) 
respectively. Lower panel: \og{q)(y) and {\og(q))(y) as computed 
using (28). The inset shows the difference (log(g)}(y) — \og{q){y) 
and demonstrates the increase in the OLR due to fiuctuations in 
q at every y. 



that it survives the mixing of PDFs from diverse 
locations into the global PDF h{r) in (51). 

As far as we are aware, bimodal RH PDFs that 
have been documented in literature are mostly derived 
from satellite data over the deep tropics (Brown and 
Zhang 1997; Zhang, Soden and Mapes 2003; Mote and 
Frey 2006; Luo, Kley and Johnson 2007). Exceptions 
are the three-dimensional and isentropic model based 
AC experiments by PBR and Yang and Pierrehumbert 
(1994) respectively: in Figure 16 in PBR and in Figure 
6 of Yang and Pierrehumbert (1994) the midlatitude 
RH PDF is clearly bimodal. 

On the other hand, the stochastic drying models 
developed by Sherwood, Kursinski and Read (2006) 
and Ryoo, Igusa and Waugh (2009) predict only 
unimodal RH PDF's. Indeed, Sherwood, Kursinski and 
Read (2006) concluded that a broad and unimodal 
RH density might be regarded as a general outcome 
of the AC model. This conclusion is inconsistent with 
our implementation of the AC model, which results in 
bimodal RH PDFs. We discuss this point further in 
section 6.3. 

Because our main focus is humidity fields on 
middleworld isentropic surfaces, we proceed with a 
data comparison by analyzing the RH on the 330K 
isentropic surface as presented in the ERA interim 
product. The data was obtained from the ECMWF 
web portal (http://data-portal.ecinwf.int/data/ 
d/interiin_daily/). We analyse daily data for the 
year 2008 given at a horizontal resolution of 1.5° x 
1.5° and at 24 vertical levels from 1000 mbar to 175 
mbar. As the RH data is only available on isobaric 
levels, we interpolate to generate RH fields on the 
330K isentrope. Further, as we are interested in the 
PDFs from the midlatitudinal portion of the isentropic 
surface, we collect data between 15° and 60° in both 
the hemispheres. 



x10 




Figure 8. The midlatitudinal RH PDF on the 330K isentrope 
as derived from the ERA interim data product for the individual 
months of 2008. Upper panel: PDFs from the two hemispheres 
plotted separately (ordering: January — back to December — 
front). Lower panel: The lines with dots are the PDFs with both 
hemispheres taken together, and the foremost bold line with star 
symbols is the mean of these monthly curves. 



A reviewer noted that the 330K isentrope 
intersects the tropopause by 60°. To alleviate concerns 
that these bimodal PDFs result from sampling distinct 
air masses separated by a mixing barrier, we have also 
constructed PDFs by considering data only from the 
troposphere, for example, from 15° and 35° in both 
hemispheres, and verified that these PDFs display the 
same qualitative form as those in Figure 8. 

As seen in Figure 8 — which shows the PDFs 
from the different months of 2008 (PDFs from the 
individual hemispheres are plotted in the upper panel, 
and a composite from both the hemispheres is shown 
in the lower panel) — there are considerable differences 
in the RH PDFs between the two hemispheres. For 
example, the dry peak in the PDFs is more prominent 
in each hemisphere's winter season. As these differences 
are not of primary interest to us we do not dwell on 
them here. There are a few values in the midlatitudes 
where the ERA dataset shows supersaturation, these 
have been set to 100% in the present calculation. 
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Naturally, this raises the secondary peak but does not 
not alter the form of the PDF. Quite clearly, the PDF 
does vary from month to month but its main features 
are robust throughout the year. More importantly, we 
observe that the ERA interim data PDF, in all months, 
is bimodal with distinct dry and moist peaks, which 
conforms with the results of the AC model in section 
5. 

6.3. Special features of the AC model in this paper 

In reconsidering the different ingredients in our 
idealized implementation of the AC model, we 
emphasize that there are at least four important 
assumptions affecting our conclusions: 

• Isentropic transport is modeled as brownian 
motion; 

• Parcels are re- moistened only at y = 0] 

• The saturation humidity is qs{y) = ^maxG""^. 

• The model neglects diffusive processes. 

The first assumption, which is equivalent to saying 
that the lagrangian velocity correlation time is zero, 
is crucial for the derivation of the Fokker-Planck 
equation, and thus for all of the results in this 
paper. Moreover, following O' Gorman and Schneider 
(2006), we emphasize that the important condensation- 
diffusion balance in (32) also relies crucially on 
the assumption of brownian motion. With nonzero 
velocity correlation time it might still possible to argue 
that mean condensation (C) is balanced by diffusive 
convergence, that is 

(C> = .,^. (52) 

But the humidity diffusivity hZq will be systematically 
less than the brownian parcel diffusivity hz^. Under- 
standing the relation between tvq and /^^ is beyond our 
scope here. 

The second assumption of a localized vapour 
source probably explains the difference between our 
conclusions and those of Sherwood, Kursinski and 
Read (2006) regarding bimodality of the RH PDFs. 
This point could be examined within the Fokker- 
Planck framework by modeling the source S in 
(10) as a re-setting, q Qsiv)^ occurring at a rate 
'?"Moist(y) throughout the domain. We do not pursue 
this generalization here. 

The third assumption, that the saturation 
humidity varies exponentially with is a consequence 
of the Clausius-Clapeyron relation and the assumption 
that the temperature gradient on a isentrope is 
uniform. We have examined alternative models of qs{y)^ 
e.g., taking qs{y) as a linear or a parabolic function of 
y while ensuring qs{y) is still monotonically decreasing 
^max at y = 0. The main qualitative features of these 
solutions shown in Figures 3, the lower panel of Figure 4 
and in Figure 6 are unchanged. The only PDF which is 
sensitive to qs{y) is the global specific humidity PDF, 
g{q) shown in the upper panel of Figure 4. Thus, as 
compared to the q~^ algebraic decay in the exponential 
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g's, a linear (parabolic) saturation qs{y) yields a uniform 
(growing) PDF as ^ ^max- These differences can be 
traced to the form of the function ys{q) in the analytic 
expression for g{q) provided by (25). 

The fourth assumption above is also crucial to 
the AC model. Haynes and Anglade (1997) estimate 
that in the upper troposphere and lower stratosphere 
the "mix-down" time for synoptic-scale fluid parcels to 
be sheared and strained to molecular-diffusive scales 
is a few weeks. This time scale is comparable to the 
mean residence time of water in the free troposphere 
(Quante and Matthias 2006). Thus molecular diffusion 
between fluid parcels may have an impact on the 
water vapour field, and in particular on its PDF. 
This important open question cannot be addressed till 
the advection-condensation model is extended to the 
advection-diffusion-condensation model. 



7. Conclusions 

We have studied the AC model forced by resetting 
the specific humidity at the southern boundary of the 
model isentrope. If the flow can be approximated by 
brownian motion (i.e., a very rapidly de-correlating 
velocity), then the PDF of the vapor field is governed 
by the Fokker-Planck equation, which can be solved 
in the limit of rapid condensation. The solution 
is illustrated for a saturation profile which decays 
exponentially along the isentrope in the meridional 
direction. Analysis of this solution, and supporting 
Monte Carlo simulations, document a non-trivial 
statistical steady state in which the southern moisture 
source is balanced by condensation throughout the 
interior of the domain. "Non-trivial" refers to the 
fact that the steady state is far from the limiting 
states of complete saturation or complete dryness, and 
instead is characterized by a spatially inhomogenous 
PDF. This non-uniformity manifests itself in the co- 
existence of very dry and nearly saturated regions 
within the domain. We find this to be encouraging as a 
strongly inhomogeneous steady state is reminiscent of 
the atmosphere. In fact, the bimodality of the global 
RH PDF predicted by the model is supported by 
midlatitudinal middleworld isentropic RH PDFs from 
the ERA interim product. 

An important property of our solution is its spatial 
structure: near the resetting source (i.e. for small y) 
the PDF is bimodal with a peak at g'min, a rapid decay 
followed by a rise to a second peak diS q ^ qs{y) ^ <7max- 
Futher from the source (i.e. for large y) the PDF 
is unimodal and terminates ai q = qs{y) <^ ^max- To 
facilitate a comparison with numerical experiments 
and observations we also derived an expression for 
the PDF as estimated over a strip of non-zero width. 
Extending the strip to cover the domain, we obtained 
the domain averaged PDF which in the case of an 
exponential saturation profile has a g^^-form. These 
features are in very good agreement with numerical 
simulations employing a large number of diffusively 
moving parcels. We also demonstrated the utility of 
the PDF in estimating mean field observables such as 
the moisture flux and condensation rate. 
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From the PDF of the specific humidity we deduced 
the PDF of the RH. Once again, the PDFs near the 
source are bimodal with a peak at the minimum RH 
for a given i.e. at rmin(^), fohowed by a decay at 
intermediate RH values and a rapid rise as r ^ 1. The 
bimo a unimodal distribution at large distances from 
the source, and here the PDF is controlled by the 6- 
function at rmin(y) which levels off to a constant as 
r ^ 1. The domain averaged RH PDF displays a mix 
of these distinct features, in particular — even though 
the global specific humidity PDF is unimodal — the 
global RH PDF is bimodal with a primary peak at 
small RH, followed by a decay at intermediate RH, and 
then a secondary peak as r ^ 1. 

Though this work deals with an idealized model, 
we believe the analysis afforded in this setting provides 
a useful understanding of advection-condensation 
and the equilibria supported by the interaction of 
these processes with a steadily maintained source 
of moisture. Indeed, with a proper identification of 
the saturation profile and the resetting protocol, our 
methodology can be applied wherever the AC model 
finds use. Also, we argue that the algorithm outlined 
here is of practical interest in that, given a saturation 
profile from a general circulation model, one can then 
predict the form of the PDF of the vapour and relative 
humidity field. In turn, this can be compared with the 
PDFs generated directly from the model itself, thus 
providing a method to validate the fidelity of moisture 
evolution in a general circulation model. 
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A. Normalization as a differential boundary 
condition 

Here we show that the normalization constraint 
obtained in the limit of rapid condensation can 
be interpreted as differential boundary condition. 
Specifically, we consider this interpretation for the 
resetting problem and pick up the story at (18) where 



S{q- 



L2 



F{q) 



y- 



(53) 



Applying the normalization (15), for G (0, L) we 
obtain 

F{q)dq=-. (54) 



y_ 

L2 



Differentiating with respect to Leibnitz's rule gives 



1 

Z2 



F{q)dq^yF[q,{y)] 

Qmin 



dy 



0. 



A second differentiation yields, 

2 



mQsiv)] 



dqs 
dy 



y 



dF fdq, 



dy V dy 



(55) 



= 0. 
(56) 



On specifying qs{y)^ we see that (56) is an ordinary 
differential equation for the function F{q). Indeed, (56) 
is the interpretation of the normalization constraint as 
a differential boundary condition. To see the utility of 
(56) in solving for F{q), consider qs{y) = e~"^, this 
yields 

dF{x) 



[2 + \og{x)]F{x) + X \og{x) - 



dx 



0. 



(57) 



and hence, F{q) = a/{Lq\og^{q)) which agrees with 
the earlier result (37). 
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